//----------------------------------*-C++-*----------------------------------//
/*!
 * \file   AngleTraits.i.hh
 * \author Jeremy Roberts
 * \date   Jul 13, 2011
 * \brief  AngleTraits inline function definitions.
 * \note   Copyright (C) 2011 Jeremy Roberts. 
 */
//---------------------------------------------------------------------------//
// $Rev::                                               $:Rev of last commit
// $Author::                                            $:Author of last commit
// $Date::                                              $:Date of last commit
//---------------------------------------------------------------------------//


#ifndef ANGLETRAITS_I_HH_
#define ANGLETRAITS_I_HH_

#include <cmath>

namespace slabtran
{

//---------------------------------------------------------------------------//
// Number of moments.

template <class D>
inline int
Moments<D>::number_moments(int N)
{
  Require(N >= 0);
  return (N + 1) * (N + 1);
}

template <>
inline int
Moments<_2D>::number_moments(int N)
{
  Require(N >= 0);
  return (N + 1) * (N + 2) / 2;
}

template <>
inline int
Moments<_1D>::number_moments(int N)
{
  Require(N >= 0);
  return N + 1;
}

//---------------------------------------------------------------------------//
// Validate moment indices.

template <class D>
inline bool
Moments<D>::valid(const int l, const int m, const int N)
{
  if ((l <= N) && (std::abs(m) <= l))
    return true;
  return false;
}

template <>
inline bool
Moments<_2D>::valid(const int l, const int m, const int N)
{
  if ((l <= N) && (std::abs(m) <= l) && ((m + l) % 2 == 0))
    return true;
  std::cout << " l=" << l << " m=" << m << std::endl;
  return false;
}

template <>
inline bool
Moments<_1D>::valid(const int l, const int m, const int N)
{
  if ((l <= N) && (m == 0))
    return true;
  return false;
}

//---------------------------------------------------------------------------//
// Generate moment cardinal index.

template <class D>
inline int
Moments<D>::index(int l, int m)
{
  Require(valid(l, m, l)); // Check m is valid for l
  // cardinal: 0       1      2      3        4      5      6      7
  //  (l,m): [(0,0)],[(1,-1),(1, 0),(1, 1)],[(2,-2),(2,-1),(2, 0),(2, 1),...
  return l * l + l + m;
}

template <>
inline int
Moments<_2D>::index(int l, int m)
{
  Require(valid(l, m, l)); // Check m is valid for l
  //
  // cardinal: 0       1      2        3      4      5        6      7
  //  (l,m): [(0,0)],[(1,-1),(1, 1)],[(2,-2),(2, 0),(2, 2)],[(3,-3),(3,-1),...
  // Only (l,m) such that l+m is even are included.  For a value of l,
  // there are l+1 terms, and the m's go as -l by 2 to l.  Before this
  // l, there are (l+0)(l+1)/2 terms.  Within the lth level, (m+l)/2 gives
  // the index, e.g. (2,2) yields 4/2 = 2.  Then (2)(3)/2=3, and 3+2=5th index.
  return l * (l + 1) / 2 + (l + m) / 2;
}

template <>
inline int
Moments<_1D>::index(int l, int m)
{
  // cardinal index is just the Legendre order for 1-d.
  return l;
}

} // end namespace slabtran

#endif /* ANGLETRAITS_I_HH_ */

//---------------------------------------------------------------------------//
//              end of AngleTraits.i.hh
//---------------------------------------------------------------------------//
